/*
 * PrimarySelector.cpp
 *
 *  Created on: 7 Sep 2011
 *      Author: Allan
 */

#include "PrimarySelector.h"

// GeoReact includes
#include "Utilities/Algorithms.h"

PrimarySelector::PrimarySelector(const ReactionSystem& reactions) :
reactions(reactions)
{}

const vector<string> PrimarySelector::SelectByFrequency() const
{
	// The species in the reactions sorted in the order of most frequent
	vector<string> sorted(reactions.GetSpecies());
	
	// Define the the comparison function
	auto compare = [&](const string& s1, const string& s2)
		{ return SpeciesFrequency(s1, reactions) > SpeciesFrequency(s2, reactions); };
	
	std::sort(sorted.begin(), sorted.end(), compare);
	
	// The number of primary species for the system of reactions
	const unsigned Nj = NumberPrimarySpecies();
	
	// The primary species as the Nj first species in (sorted)
	vector<string> primarySpecies(sorted.begin(), sorted.begin() + Nj);
	
	return primarySpecies;
}

const vector<string> PrimarySelector::SelectByConcentration(const VectorXd& n) const
{
	// The species in the reactions 
	const vector<string>& species = reactions.GetSpecies();
	
	// The indexes of the above species
	const vector<Index>& indexes = reactions.GetIndexes();
	
	// The composition of the above species
	VectorXd ni(indexes.size());
	
	for(unsigned i = 0; i < indexes.size(); ++i)
		ni[i] = n[indexes[i]];
	
	// The species in the reactions sorted in the order of most concentrated
	vector<string> sorted(species);

	// Define the the comparison function
	auto compare = [&](const string& s1, const string& s2)
		{ return ni[IndexOf(s1, species)] > ni[IndexOf(s2, species)]; };
	
	std::sort(sorted.begin(), sorted.end(), compare);
	
	// The number of primary species for the system of reactions
	const unsigned Nj = NumberPrimarySpecies();
	
	// The primary species as the Nj first species in (sorted)
	vector<string> primarySpecies(sorted.begin(), sorted.begin() + Nj);
	
	return primarySpecies;
}

const unsigned PrimarySelector::NumberPrimarySpecies() const
{
	return reactions.GetSpecies().size() - reactions.GetReactions().size();
}

const unsigned PrimarySelector::SpeciesFrequency(const string& species, const ReactionSystem& reactions) const
{
	unsigned count = 0;
	
	BOOST_FOREACH(const Reaction& reaction, reactions.GetReactions())
		if(IsContained(species, reaction.GetSpecies()))
			++count;

	return count;
}
